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ABSTRACT: We use Quantum Monte Carlo to evaluate the conductivity 
a of the 2-dimensional disordered boson Hubbard model at the superfluid-bose 
glass phase boundary At the critical point for particle density p = 0.5, we find 
cTc = (0.45i 0.07)crQ, where aq = e1/h from a finite size scaling analysis of the 
superfiuid density We obtain ac = (0.47— 0.08)(tq from a direct calculation of the 
current-current correlation function. Simulations at the critical points for other 
particle densities, p = 0.75 and 1.0, give similar values for a. We discuss possible 
origins of the difference in this value from that recently obtained by other numerical 
approaches. 
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Introduction 



The interplay between interactions and disorder has been a compelling field of research 
over the last decade. Particular interest has focused on two dimensions, both because ran- 
domness alone marginally localizes the eigenstates, and also because of a set of fascinating 
experiments on the superconducting-insulator transition. Early studies"*^ of the sheet resis- 
tance as a function of temperature for granular systems of Sn and Ga on glass substrates 
suggested that the resistance was the only important variable controlling whether the sys- 
tem went insulating or superconducting at low temperatures. Films differing significantly 
in thickness or other characteristics, but having similar resistance, would end up in the 
same low T ground state phase. Furthermore, it was found that the transition to the super- 
conducting state consistently took place at a resistance close to h/^e^ = BABkQ. On the 
other hand, experiments^ on homogeneously disordered Bi and Pb films on Ge substrates 
explicitly employed the thickness as a control parameter- for thin films the resistivity is 
characteristic of an insulator, i.e. p increases as the temperature is lowered, while for thick 
films a superconducting phase transition is observed. In between is a separatrix which has 
the property that, roughly speaking, p neither diverges nor goes to zero as T ^ 0. For 
Bi the value of the resistance at the critical thickness is similar to that found for granular 
systems, namely close to 6A5K^. For Pb the value is somewhat larger, around 9.5K^l. 
Since these materials were homogeneously disordered, it was suggested that the universal 
resistance phenomenon must be attributed to a more general principle than those originally 
suggested to explain the behavior of strictly granular materials. 

Meanwhile the mechanism for the destruction of superconductivity was also explored. 
In granular Al, Pb, and Sn films, it was found^ that the superconducting gap and Tc varied 
rather little with thickness. The insulating transition was driven by long resistive tails 
which gradually broadened until they reached the size of the gap itself. For homogeneous 
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Pb films, on the other hand, and the gap were driven to zero by the disorder. This work 
emphasized the possibifity of destroying superconductivity either through a reduction of 
the amphtude of the pair wave function or alternately through the destruction of phase 
coherence.^ 

Since then, many experimental groups have explored related effects for different ma- 
terials and different control parameters. Ga was observed^ to have a critical resistance of 
QKVt. In MoC films, again using the thickness to dial through the transition,^ a thresh- 
old resistance in the range 2.8 — 3.5X1) was measured. In doped, amorphous semiconductor 
films (InOa;) of fixed thickness, the degree of microscopic disorder has been used to dial 
the transition.^ Also in InOx (i?c ~ hKVL) a magnetic field can drive^ the systems across 
the superconducting-insulator phase boundary. Finally, in Josephson junction arrays^ the 
ratio of charging to Josephson energies is similarly tuned. 

While the experiments all describe qualitatively similar behavior, the quantitative sit- 
uation is still developing. Not only does the "universal conductance" in fact vary, but 
there can be non-trivial structure in the curves near the separatrix as the temperature 
is lowered. In particular, a "re-entrant" phenomenon is observed in which the resistivity 
dips as if the film were about to go superconducting, but then the transition is usurped 
by the formation of an insulating state and the resistance rises as T is decreased further. 
Further questions concern whether the experiments are really in the critical regime or not. 
That the experimentally accessible temperature ranges, T > 0.5K, may not be sufficiently 
low is an issue that has been raised by, among other things, the existence of structure in 
the curves as T is reduced. This is a crucial question in determining the appropriate model 
for the critical phenomena, as we shall discuss further. 

There have been a number of theoretical efforts to understand these phenomena. -"^^ A 
particularly interesting set of ideas has centered on the proposal that, despite the fact that 
the underlying degrees of freedom are fermionic, the appropriate model is one of disordered, 
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interacting bosons/^'^^ 

{ij) i i 

Here a]^ is a boson creation operator at site i, ni — a^ai is the number operator, F is a 
soft-core on-site interaction, and is a random local chemical potential, which we chose to 
be uniformly distributed between —A and A. The physical motivation for a boson model 
is provided by imagining one has pre-formed Cooper pairs above the superconducting 
transition temperature, i.e. the fermions first condense into a set of interacting bosons 
that lack phase coherence. At lower temperatures phase coherence is established, and the 
system becomes superconducting. Such a view is particularly natural for granular systems 
where one can imagine Cooper pairs forming on individual grains. Cha et.al}'^ have 
presented a more careful characterization of this general scenario by comparing correlation 
lengths for the various operator expectation values associated with the propagation of two 
electrons both together and independently. The fundamental physics is that on the scale 
of the diverging pair-correlation length Cooper pairs appear as point bosons. 

However, the justification for considering the Hamiltonian Eq. 1 is by no means com- 
pletely qualitative. The renormalization group calculation of Giamarchi and Shultz^^ pro- 
vides at least one explicit theoretical demonstration that the universality classes describing 
the superconducting transition of fermions with an attractive interaction and the superfluid 
transition for disordered bosons are identical. More detailed ea;^?enmenta/ justification of 
a picture of pre-formed bosons has been provided by the work of Paalanen et.alM, which 
lends rather compelling support to a picture where the amplitude of the Cooper pair wave 
function is finite on both sides of Tc, and the superconducting transition is indeed de- 
scribed in terms of a loss of phase coherence between pairs rather than the breaking of the 
pairs themselves. Even so, it must be pointed out that tunneling experiments^' on ho- 
mogeneously disordered Pb films have indicated a vanishing of the gap at Tc- It therefore 
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may be that granular and Josephson junction systems are more appropriate realizations 
of bosonic models. 

Here we will describe the results of quantum simulations combined with a finite size 
scaling analysis to determine the conductivity at the superfluid-insulator transition in the 
boson-Hubbard model. The organization of this paper is as follows: We will first review 
briefly some of the analytic and numerical work to put our calculation in context. Next we 
describe our monte carlo method and discuss how the conductivity can be obtained either 
from a Kubo formula for the current-current correlation function, or, alternately, from the 
rigidity of the system to phase changes. We then detail our numerical results, beginning 
with the determination of the critical point for the bose glass-superfluid transition, and 
continuing with the analysis yielding the conductivity. The two techniques yield values 
which are the same to within our estimated error bars. However, this number differs 
substantially from that obtained previously. -"^^'^^ We conclude with a discussion of some 
possible sources of this discrepancy. 

Fisher et.alM have qualitatively mapped out the ground state phase diagram of the 
boson-Hubbard model. In the clean limit, a gapless superfluid phase exists for all non- 
integer fillings. At commensurate densities, however, the interactions freeze the bosons 
into a gapped Mott insulating (MI) phase for sufliciently strong coupling. Increasing 
hybridization t will eventually wash out the gap and drive the system from insulator to 
superfluid. It was argued that this "coupling-driven" transition is in the universality class 
of the classical d+1 dimensional XY model. By contrast, changing the density away from 
an integer number of bosons per site can also induce superfluidity, but here the transition 
is mean fleld in character. When disorder is added, it was suggested ■'^■'^ that a third, "Bose 
glass," phase appears. This phase is characterized by the absence of a gap, but also by a 
vanishing superfluid density. While all non-integer boson densities were superfluid prior 
to the introduction of randomness, an incommensurate insulating phase is now possible 
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as the disorder increases. This work^^'^'^ also predicted values or bounds for the critical 
exponents, some of which have been verified experimentally.^ 

The ground state phase diagram was subsequently mapped out numerically^^'^^''^^'^^''^'^ 
in one and two dimensions and also studied by the Bethe Ansatz.^'^ Quantitative values for 
the coupling required to localize the bosons into the MI phase, were determined, ^^'^^'^'^ and 
the prediction of mean-field exponents for the density controlled transition was verified. -"^^ 
In the presence of disorder, the basic picture of the formation of a bose glass phase was sub- 
stantiated, although a variety of unexpected re-entrant phenomena were also oberved.^^'^^ 

Numerical and analytic efforts have more recently turned to the transport properties 
for two dimensional lattices, and, in particular, the evaluation of the conductivity. Runge 
used exact diagonalization techniques, combined with finite size scaling to extract the 
conductivity of the disordered, hard-core {U — oo) model. He found Uc = (0.15^ 0.01)a"Q. 
In a set of papers, ■'^^'^'''^^ the 3-d XY model and its Villain variant were studied with and 
without disorder and with short and long range potentials. It was argued that this model is 
in the same universality class as the disordered boson-Hubbard Hamiltonian. Cha et.al.^^ 
obtained the conductivity first in the clean system with short-range interactions. The 
value (Tc = (0.29i 0.02)(Tg found by Monte Carlo compared favorably with an analysis 
based on a 1/N expansion which gave Uc = 0.251(Tg. Sorensen et.al}'^ then showed that 
the addition of randomness results in a 5ma//er value ctc = (0.14^ 0.01)crQ, that is, farther 
from the experimentally realized numbers. Finally,^^ the disordered model with a long 
range Coulomb potential was studied. Including these interactions was found to push the 
conductivity back up to ctc = (0.55^ 0.01)crQ. This final value is certainly within a factor 
of two or so of the experiments, and possibly substantially closer especially considering 
uncertainties associated with the precise low temperature experimental values. 

None of these studies were of the Hamiltonian Eq. 1. It is of interest to compute the 
properties of the boson-Hubbard model directly, including the effect of the number fiuctu- 
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ations ignored in the mapping to the XY model and also in the hard-core diagonalization 
methods. One motivation is to test the arguments suggesting the universality classes are 
identical. In addition, if the experiments are not in the critical regime, then non-universal 
quantities become of interest, and the predictions of the original model are essential. 

Monte Carlo and Finite Size Scaling Methods 

Here we present a brief discussion of our numerical approach. More detailed descrip- 
tions have recently appeared. ^^'^-"^ We begin by expressing the partition function as a path 
integral. In order to do this, we discretize the imaginary time (3 — L^-Ar and use the 
Trotter approximation^^ to decompose the imaginary time evolution operator. 

Z = Tr e-^^ = Tr [e"^^^]^- ~ 1^ [JI (2) 

i 

This is a well-controlled procedure since one can explicitly calculate at different At and 
use well understood techniques'^ to extrapolate to At = 0. We now insert complete sets 
of states to express Z as a sum over a classical occupation number field n{l, t). 

Z= J2 {n{f,lK^^'\n{t,2){n{r,2\e-^^^\n{t,3) .... (3) 

{n{l,r)} 

We have thus written our d dimensional quantum mechanical trace as a classical statistical 
mechanics problem in d + 1 dimensions. The classical degrees of freedom to be sampled 
are the space-imaginary time dependent boson occupation number field. Due to particle 
number conservation in H, the allowed configurations of this field trace out "world-lines" 
as the particles propagate in r. The "Boltzmann weight" of a particular configuration 
is the product of matrix elements, and can be sampled with standard stochastic tech- 
niques. Various choices are possible for the decomposition of H into Hj^. We chose the 
"checkerboard" break-up in these studies.'^ 

Expectation values are constructed in a similar way. Of particular interest to us will 

7 



(4) 



be the true (paramagnetic) current-current correlation function Jxxi'^) 

Jix{T) = {jUr)fxm 

fx{r) = e^-me-^' 
i 

and its fourier transform 

Jix{iwm)= f dT{fx{T)jlmei^-\ (5) 

Jo 

If we were able to evaluate these quantities at real frequencies, a;, we could then 
compute the frequency dependent conductivity from the Kubo formula 30,31,32 

a{uj) = <jQ — [-jUu:)-{kx)]. (6) 

UJ 

Here kx is the kinetic energy in one of the two equivalent lattice directions. The dc 
conductivity is then obtained from aj^, = (j{u 0). However, it is well-known that 
the analytic continuation of quantities like the conductivity to real frequencies is a subtle 
problem. "^"^ Following Sorensen et.al.}'^ our approach will be to exhibit that for small 
Matsubara frequencies {uim — 2rmr/(3), the conductivity obeys a{iu3m) — + |<^m|T'c)- 

With the empirical observation that our data fit this analytic form, the continuation to 
(7{uj+i6) = a>K/(l— za;Tc) is then possible. A direct calculation of cr, avoiding the assumption 
of the Drude form, similar to the computation of real frequency quantities in fermion 
simulations,'^^ would be useful. 

There is an alternate approach to obtaining the conductivity based, instead, on the 
response of the system to a phase twist. We define a "pseudocurrent" operator and its 
correlation function by 

^(^) = ^J20xir + r')jx{r')+jy{r + r')jy{r')), 

N, (7) 
jx{r) = r + 1) - x{i, r)]. 

i=l 
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Here x{i, r) is the x-component of the position of boson i at time shce r, A^^ is the total 
number of bosons and L the extent of the lattice in each spatial direction. Therefore, jxiT) 
measures the number of bosons moving in the positive x direction minus those moving in 
the negative direction at time slice r, i.e. the net boson flux per time step at time slice r 
in the x direction. Similar definitions apply for the pseudocurrent in the y direction. This 
pseudocurrent operator was originally introduced-*^^ to measure the mean square winding 
number and from that the superfluid density. 

T 

J{u; ^ 0) = (8) 

Defining ps(<^) = i7(<^)/2t, we can show that the frequency dependent conductivity can be 
written as 

^ = 4nt^. (9) 

Again a Drude form must be assumed here to carry out the analytic continuation. 

Clearly, these two approaches are related, since they both yield the conductivity. In- 
deed, a detailed discussion of the connections between the response of the free energy to 
changes in the boundary conditions and the current-current correlation function has re- 
cently appeared."^^ However, despite these general relationships, the formulations of Eqs. 4- 
6 and 7-9 allow us to evaluate a in two quite different ways. In particular, the current- 
current correlation function is a simple operator expectation value for the model described 
by the boson-Hubbard Hamiltonian. Meanwhile, the pseudocurrent analysis is based on 
topological properties of the boson world lines in our path integral representation of the 
partition function. Thus we regard the two measurements as rather independent determi- 
nations of the conductivity and a check for the self consistency of our simulations. 

It is useful to contrast our Monte Carlo approach with the two other numerical tech- 
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niques already used in evaluating the conductivity. As discussed above, Runge^" used an 
exact diagonalization method for the U = oo model. The advantage of that approach is, 
most importantly, that real time correlation functions can be directly inferred, without the 
need for any analytic continuation. A second advantage is that there are no statistical er- 
rors from Monte Carlo sampling, only fluctuations associated with the disorder averaging. 
On the other hand, the technique is limited to rather small lattices (2 x 2 up to 5 x 5) 
since the Hilbert space dimension grows exponentially with the number of sites. Indeed, 
the Hilbert space also grows rapidly with the number of bosons allowed per site, so that 
in practice it is necessary to consider the hard core case where the site occupations are 
limited to 0, 1. While this prevents the study of the transition as a function of interaction 
strength, as we will do, one is able to evaluate the conductivity for C/ = oo at a critical 
point accessed by changing the density. 

Meanwhile, Monte Carlo simulations have also been conducted of a spin model (the 
classical 3-d XY model) argued to be in the same universality class as the boson Hub- 
bard Hamiltonian. These studies are very close in spirit to the ones described here, since 
they are also simulations of a classical model in one higher dimension than the original 
quantum Hamiltonian. Unlike diagonalization, they are subject to error bars associated 
with statistical sampling and have to argue the analytic continuation, but they are able to 
study lattices of order 100 sites. There are two differences with our approach. The action 
in our simulation, the product of matrix elements of Eq. 3 is rigorously appropriate to the 
boson-Hubbard Hamiltonian. For example, the mapping to the spin-1/2 model neglects 
boson number fluctuations. Presumably, the advantage of the spin model approach is that 
the action is somewhat simpler than that arising from the quantum^classical "world-line" 
description detailed above. On the other hand, there are evidently assumptions associated 
with the universality class of the transition and the relevant degrees of freedom. 
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Numerical Results 



In order to evaluate cr, we must first determine the superfluid-insulator critical point. 
Here we closely follow the finite size scaling procedure of Ref. 

According to two parameter finite size scaling, physical quantities, such as ps, on 
different size lattices, of linear extent L, satisfy^^ 

Ps = L^f{aLh,PL-'). (10) 

Here a = 2 — d — (3 is the inverse temperature, and 5 = {V — Vc) jVc measures the 
distance to the critical point. The function / is universal but the metric factor a is not. 
The dynamic critical exponent z is predicted"^^ to have the value 2; = 2, and our system 
is two dimensional giving a = —2. Notice that there is no nonuniversal metric factor 
for the second argument of the function /. Therefore, by keeping the second argument, 
PL~^, fixed and plotting L'^ps versus V for different lattice sizes, all the curves should 
intersect at the critical value Vq. We chose the inverse temperature f3 = (1/4)L^At for 
different lattices L, i.e. the same aspect ratio as in Ref.l7. Data for L^s, obtained from 
the pseudocurrent correlation function via the procedure described in Eqs. 7-9, is shown 
in Fig. 1. We see that the curves from different lattice sizes converge at a value V/t !v 7.0, 
but then do not fan out. The reason for our inability to resolve the exact point at which 
these curves intersect is the uncertainty introduced by the extrapolation of the superfiuid 
density, to a; = 0. Because of this, we use Fig. 1 only to give us a rough idea of where 
the critical point is. Note that at the critical point, plots of p{ijj)/uj versus uj for different 
lattice sizes should collapse on a single curve (see Eq. 9). We use this collapse of the data 
onto a single curve to determine the critical point more accurately than possible with figure 
1. 

In our simulations at half filling, p = 0.5, the strength of the disorder was kept constant 
at A = 6, and the coupling, V, varied to find the critical point. The number of disorder 
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realizations we for most of our simulations is of the order of 100 (see the figures for 
more details). For the 8x8 and 10 x 10 lattices we did about 40000 thermalization 
and 100000 measurement sweeps, for the 12 x 12 lattice we did 100000 thermalization and 
20000 measurement sweeps. We found that the large number of sweeps was necessary for 
thermalization and good statistics. 

Figures 2a-c show plots of p{ijj)/uj versus uj for different values of the coupling, V/t = 
6.5,7.0,8.0. We see that the best data collapse is at V/t — 7.0. Furthermore, the plots 
show that the values 6.5 and 8.0 actually bracket the critical region. This is because for 
V/t = 6.5 the data for 8 x 8 x 16 lie slightly below those for 10 x 10 x 25, as one would 
expect if the system is in the superconducting phase, while the pattern is reversed for 
V/t = 8.0, as should happen when the system goes into the insulating phase. To obtain 
the conductivity, we fit the Monte Carlo data using (t{u)) = a/(l + biVm) (solid lines in the 
figures). If we use the data for V/t = 6.5 we find a ^ 0.5crg, whereas for V/t = 7.4 (not 
shown) we get a = O.SScrg. The best data collapse, for V/t = 7.0 gives a = 0.4crQ. To 
check the effect of the finite time step errors, which are known to be O(At^), we redid the 
simulation for V/t = 7.0, with the same /? but twice as many time slices, i.e. half the At 
as before. The results are shown in Fig. 3 and give a = 0.44crg. Extrapolating to At — > 
gives a/aq = 0.47± 0.08. 

We now compute the conductivity from the true current-current correlation function 
using the Kubo relation by taking the limit 17(0; 0) in Eq. 9. In Fig. 4 we show the 
fourier transform of the current-current correlation function, J^iuj), as a function of uj at 
the critical point V = 7.0 for A = 6. for lattice sizes iV = 8 x 8, 10 x 10, 12 x 12. To obtain 
the conductivity, we fit the data with the Drude form 

Jiu) - Jiu = 0) = (11) 

According to the Kubo formula, Eq. 6, the kinetic energy should be subtracted from 
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the frequency dependent current-current correlation function to obtain the conductivity. 
Scalapino et.al.^^'^'^ demonstrated recently that whether or not the q = current-current 
correlations extrapolate to the kinetic energy A; as a; — > can be used as a signature of 
the superfluid-insulator transition in the fermion Hubbard model. We have verified that 
this is true for bosons too because our independent measurements of the current and the 
kinetic energy do correctly signal this transition in agreement with other order parameters 
such as the superfluid density and gap. Indeed, this analysis provides a separate check on 
our code. Fig. 5 is the same as Fig. 4 but with twice as many time slices. 

Table 1 shows the values for the fitting parameters A and T, as well as the associated 
conductivity 2itA = 0.38-0.46 with the range in values representing changes due to different 
spatial and imaginary time lattice geometries. These numbers are in agreement to within 
our systematic and statistical uncertainties with those found from the independent analysis 
of the pseudocurrent in the previous section. 

To check universality, we did similar simulations for p — 0.75 and p — 1- The con- 
ductivity for p = 0.75 near the transition is shown in Fig. 6. These results were obtained 
using simulations of 5 — 20 thousand sweeps for thermalization and 20 — 40 thousand for 
averaging. Fifty realizations of the 8x8 lattice and fifteen realizations of the 10 x 10 lattice 
were used. In these simulations tA/V was held fixed at 3/7 while varying A and V/t. The 
transition was determined by searching for where the conductivity changes from increasing 
with size to decreasing with size, i.e. as in Fig. 2a-c. In the Bose glass phase near the tran- 
sition we find that for the small lattices we are using it is difficult to distinguish between 
data which collapse on the same curve and data where the conductivity decreases with 
size. This may be one of the reasons why previous studies have found lower values for a. 
Fits of the data in Fig. 6 to the form cr(a;) = a/(l -|- bujm) show that cr/aq varies from 0.36 
to 0.38. The range of values depends on the number of data points used in the fit. Using 
all the data gives the lower values; using only the two points closest to ojm = (which have 
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larger statistical error bars) gives the higher values. At p = 1 we obtain similar results as 
shown in Fig. 7. This is an important fact since it is possible that the Bose glass phase 
would not completely cover the Mott lobe, and thus at integer filling one might find a 
phase transition from the superconductor directly into the Mott insulator. The data were 
obtained with 10, 000 sweeps for thermalization and 40, 000 sweeps for averaging. From 
60 to 170 realizations were used for lattice sizes of L = 8 and L = 10 and fixed V/t = 10. 
These results show reasonable scaling. Again, estimates for a/aq are in the range 0.3 to 
0.4. The conductivity results for p — 0.75, 1 are consistent with those for p — 0.5. 

Our most accurate data are for half filling, and these were done on the Connection 
Machine CM5 with 32 and 64 processing nodes. The simulations for p = 0.75, 1 were done 
on Silicon Graphics workstations using a program written independently from the one used 
for half-filling. 

Conclusions 

In this paper we have measured the conductivity of the disordered boson-Hubbard 
model. We obtain a value cr = 0.46^ O.OSug which differs from those previously reported 
for studies of a classical spin model^^ and diagonalization of hard-core boson systems on 
small lattices. One possible source of this discrepancy is that we noticed in our data that 
the scaling curves which are extrapolated to zero frequency to obtain the conductivity, are 
still curving upwards towards large values for the larger lattices. It is therefore possible 
that for larger lattices than examined here (12 x 12) the conductivity could be even larger. 

The precise source of this discrepancy is still under investigation. We believe that 
by evaluating the conductivity by two rather different approaches we have reduced the 
possibility of trivial questions of normalization. This still leaves the possibility of various 
numerical uncertainties, including questions concerning whether we average over a suffi- 
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cient number of disorder realizations, the nature of Trotter errors, etc. We have exphcitly 
checked this latter point by running at two values of At and comparing the results. We 
found that this source of systematic error is small. We also carefully checked the equili- 
bration of our lattice. The biggest source of uncertainty in our calculation is clearly in our 
determination of the critical point. From Fig. 2, it seems clear that our value for cr is rather 
sensitive to this choice, and a range of values is possible. However, we do not believe that 
the previously reported numbers are consistent with our data. This would require a choice 
of Vc clearly inconsistent with any sort of scaling plot. We are in the process of carrying 
out some exact diagonalization studies to pin down the source of the disagreement with 
QMC calculations. 
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FIGURE CAPTIONS 



1. The scaling variable {L Ps) vs. the interaction strength, V. The superfluid density, ps 
is obtained from the pseudocurrent correlation function via the procedure described in 
Eqs. 7-8. The convergence region of the different curves gives the approximate critical 
value for V. 

2. Plots of ps{u))/u) versus ui for different values of the coupling, V/t = 6.5, 7.0, 8.0. The 
best data collapse indicates the critical point, and that happens for V/t = 7.0. 

3. Plot of ps{oj)/oJ versus uj for coupling V/t — 7.0. Here the number of time steps has 
been doubled to check for Trotter errors. 

4. Jioj) as a function of uj at the critical point V = 3.5 for A = 6 and lattice sizes 
A?" = 8 X 8, 10 X 10, 12 X 12. The good data collapse again indicates that the system is 
at the critical point. 

5. J{ijj) as a function of u at the critical point V = 3.5 for A = 6, for lattice sizes 
iV = 8 X 8, 10 X 10, 12 X 12. Here the number of time steps has been doubled (at fixed 
P) to check for Trotter errors. 

6. Scaling plot as in figure 2 but at three quarter filling, cr/ug ~ 0.37. 

7. Same as figure 6 but at p = 1. a/aq ^ 0.35. 
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Table I. 





Ar 


.4 


T 


a 


8x8 


0.222 


0.061 


0.33 


0.38 


10 X 10 


0.222 


0.065 


0.38 


0.41 


12 X 12 


0.222 


0.068 


0.40 


0.43 


8x8 


0.111 


0.067 


0.36 


0.42 


10 X 10 


0.111 


0.072 


0.42 


0.45 



Table 1 Caption: The values of the fitting parameters A and T used to obtain the smooth 
curves in Figs. 4a-c, and the associated number for the conductivity. A and T were 
obtained by fitting the analytic form Eq. 11 to the Monte Carlo data at a; = 1 and 
a; = 2, a process which results in an excellent fit even out to much higher frequencies. 
The imaginary time dimension has been chosen as L = N'^/A for At = 2/9 and chosen as 
L = for At = 1/9, so that the inverse temperature /? is the same for the two different 

discretization intervals. We have also done least squares fits to the entire frequency range, 
a process which generally gave values for the conductivity which were higher by 10% or 
so. However, we believe that fitting to the low frequency data is the more appropriate 
procedure, since the analytic form is expected to be most valid there. 
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